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Abstract 

Simple models of magnetic field generation by convection in rotating spherical shells exhibit properties 
resembling those observed on the sun. The assumption of the Boussinesq approximation made in these 
jnodels prevents a realistic description of the solar cycle, but through a physically motivated change in the 
boundary condition for the differential rotation the propagation of dynamo waves towards higher latitudes 
pan be reversed at least at low latitudes. 
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.1. Introduction 

A well known difficulty in modeling the solar dynamo is the 
'fact that the dynamo waves which describe the nearly time pe- 
riodic dynamics of the magnetic field tend to propagate from 
lower to higher latitudes instead in the opposite sense as ob- 
served on the sun. This effect is well known in mean field mod- 
pis of the solar cycle, see for example Stix (1976), but it is also 
observed in numerical solutions of convection-driven dynamos 
in rotating spherical fluid shells which are supposed to model 
processes in the solar convection zone. This shortcoming and 
jothers are caused by the inadequate representation of the com- 
pressibility of the solar atmosphere. Even with the huge power 
pf modern supercomputers it is not yet possible to resolve ade- 
quately the dynamics of convection in the presence of the large 
'density variation between bottom and top of the solar convec- 
tion zone and to model appropriately the dependence of density 
pn pressure. 

' In their early direct numerical models of solar convection and 
magnetic field generation. Oilman and Miller (1981) assumed 
the Boussinesq approximation in which the fluid is regarded 
as incompressible except in connection with the gravity term 
where the temperature dependence of the density is taken into 
account. This approximation eliminates the need for a sepa- 
rate equation of state and leads to a system of equations de- 
scribing long period processes while the short period acoustic 
modes no longer enter the analysis. The same effect is obtained 
in the anelastic approximation in which the horizontally aver- 
aged density variation is taken into account, but the fluctuating 
component of the density is still only represented in the gravity 
term. For applications of the anelastic approximation in models 



of the solar convection zone see Oilman and Olatzmaier (1981) 
and ElUot et al. (2000). 

Convection in rotating spherical fluid shells heated from be- 
low is always associated with a differential rotation generated 
by the Reynolds stresses of convection. An analytical model 
demonstrating the preference of banana shaped convection cells 
girdling the equator and the associated solar like differential ro- 
tation was presented by Busse (1970, 1973). While the ana- 
lytical solution for stress-free boundaries exhibits a depth inde- 
pendent differential rotation, a differential rotation decreasing 
with depth is always found in fully nonlinear numerical mod- 
els. This property together with the fact that the differential 
rotation reaches its maximum at the equator is responsible for 
the propagation of the dynamo waves towards higher latitudes 
(Yoshimura, 1975). 

The solar differential rotation also decreases with depth 
throughout most of the convection zone with the exception of 
the tachocline region near its bottom and a region near its sur- 
face as indicated in figure 1(b). There is no general agreement 
about the origin of the upper 30Mm deep layer in which the 
differential rotation increases with depth. Here we assume that 
it is caused by supergranular convection that is characterized 
by a strong asymmetry between rising hot plasma and descend- 
ing cool plasma. This type of convection has been modeled by 
hexagonal convection cells in the presence of rotation (Busse, 
2007). As shown in this paper the asymmetry between rising 
and descending flow in hexagonal convection cells does indeed 
generate a differential rotation that increases with depth. We use 
this dynamical property of convection as motivation to modify 
the usually assumed stress-free boundary condition. Solely for 
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the differential rotation we apply the condition given in expres- 
sion (All) of the Appendix. The resulting profiles shown in 
the example of figure 1(a) are still not very solar like, but they 
show an increase with depth of the differential rotation in the 
upper layer and tend to generate a more solar like behavior in 
the numerical simulations, as discussed below. 

2. Simple Convection-Driven Dynamos 

Because convection driven dynamo solutions strongly de- 
pend on the Prandtl number P and on the magnetic Prandtl num- 
ber P,„ for values of the order unity and because even bistabil- 
ity has been found in this parameter range (Simitev and Busse, 
2009) we use a numerical model with a minimum number of 
external parameters as outlined in the Appendix. In this way it 
is possible to cover a relevant parameter region in sufficient de- 
tail. In addition to the two Prandtl numbers and the parameter 
jS introduced in the boundary condition (All), only the rotation 
parameter t and the Rayleigh number R must be considered. 

In choosing parameter values for a numerical solar model the 
task is most easily accomplished in the case of t since only a 
value of the (eddy-) viscosity must be chosen. Using a com- 
monly accepted eddy viscosity of the order IQ^m^/s (Oilman, 
1983) we findr ^ 2 - 10^ It is more difficult to select an appro- 
priate range for the Rayleigh number. We shall consider values 
that exceed the critical values for the onset of dynamo action by 
less than a factor of two. Otherwise convection motions become 
too chaotic and the structures of the dynamo solutions become 
less regular and more difficult to understand. 

In figure 2 the properties of a convection-driven dynamo 
are indicated by following the nearly periodic magnetic field 
through half a magnetic cycle. The left column shows on the 
left side lines of constant magnetic flux density where the 
bar indicates the average over the azimuthal coordinate (p. On 




Figure 3: (color online) Contour lines of the radial velocity at r = r, + 0.5 
and the poloidal streamlines rdv/dtp in the equatorial plane in the case P = 1 , 
r = 2000, fi = 120000, P,„ = 4, ^8 = 1.5. 

the right side of the meridional cross-section lines of the ax- 
isymmertic poloidal field are shown. The next column indi- 
cates the places where the toroidal field reaches its maximum 
strength near the surface. Since the toroidal field provides by 
far the major contribution to the horizontal magnetic flux den- 
sity, eruptions of magnetic bipolar regions occur most likely 
at those places in the case of the sun. The modification of 
the stress-free boundary condition given by expression (All) 
causes the differential rotation to increase with depth in the up- 
permost layer of the shell. This is turn causes the "magnetic 
bipolar regions" shown in the second column of figure 2 to 
drift in a solar-like fashion from higher latitudes towards the 
equator in agreement with solar observations. It is also remark- 
able that non-axisymmetric components characterized by the 
azimuthal wavenumber m = 1 of the magnetic field dominate 
over the axisymmetric component as is evident from the third 
column which shows the strength of the radial magnetic field 
at the surface of the sphere. The superposition of the axis- and 
non-axisymmetric components actually represent the property 
that the dynamo process is nearly confined to one meridional 
hemisphere while the other hemisphere is almost field free. 
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Figure 2: (color online)Approximately half a period of a dynamo oscillation in the case t = 2000, R = 120000, P = 1, P,„ = 4, /? = 1.5. The first column 
shows meridional lines of constant on the left and poloidal fieldlines, r sinddh/dd on the right. The second column shows fines of constant dgjdd at r = 0.9 
corresponding to -0.9, -0.8, -0.7, 0.7, 0.8, 0.9 of tfie maximum absolute value of dg/d8, and the third column shows lines of constant 6^ at r = r^. Tfie last column 
sfiows Re{dg"'~^ / dd) on the left and lm(dg"'~^ / dd) on the right. The five rows are separated equidistantly in time by Af = 0.0224. 



The solar evidence for a non-axisymmetric magnetic field is 
provided by the phenomenon of active longitudes (see Usoskin 
et al. (2007) and references therein). As in the case of those lon- 
gitudes, the field shown in the third column of figure 2 hardly 
propagates throughout the time frame of the figure. The angu- 
lar velocity of the prograde propagation of the pattern is less 
than 2n which is similar to the propagation of the convection 
pattern. In contrast to the cyclical nature of the magnetic field 
the convection pattern is nearly steady as indicated in figure 3. 
Individual convection columns may sometimes breakup in that 
the outer part tends to propagate a bit faster than the inner part 



which essentially remains steady in the rotating frame. 

In figure 4 half a magnetic cycle is presented for a higher 
Rayleigh number than in the case of figure 2. As is evident 
from the 3 column of figure 4 the structure of the magnetic field 
is not antisymmetric with respect to the equatorial plane. In the 
southern hemisphere the m = 2 -component of the magnetic 
field dominates in this particular cycle. At other times patterns 
that are more antisymmetric and similar to that shown in figure 
2 are found. The amplitude of the magnetic field varies strongly 
from cycle to cycle as shown in figure 5 while amplitude of 
convection exhibits rather small variations. 
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Figure 4: (color online)Nearly half a period of a dynamo oscillation in the case r = 2000, R = 140000, P = 1, P,„ = 3.5, /3 = 1.0. The first column shows 
meridional lines of constant on the left and poloidal fieldhnes, rsmddh/d0 on the right. The second column shows contour lines of |6;,oricOTiin;|sgn(B^) at r = 0.9 
corresponding to -0.9, -0.8, -0.7, 0.7, 0.8, 0.9 of the maximum value of |B;,ori;o/im/l, and the third column shows lines of constant at r = To. The last column shows 
Re(dg"'~^ / d6) on the left and lm{dg"'~^ / 06) on the right. The four rows ai'e separated equidistantly in time by Af = 0.0168. 



It is tempting to present visualizations similar to the famous 
butterfly diagrams of solar cycles that have recorded observa- 
tions of sunspots for many decades. Since we associate the 
occurrence of sunspots with the strength of the horizontal com- 
ponent of the magnetic field in the sun near the surface, the am- 
plitudes of the maximum values of the horizontal components 
could be plotted as function of latitude and time. But this pro- 
cedure could not be realized due to limitations of data storage. 
Instead we have restricted the attention to the main contributer 
to the azimuthal component of the magnetic field. Thus, in the 
top panel of figure 6 the extremal azimuthal component of the 
magnetic field corresponding to the wave numbers m - and 
m = l,i.e. Z?^ -I- |Z?™"'|sgn(Z?^), is plotted as function of time and 
latitude. In contrast to solar butterfly diagrams the movements 
to higher latitudes still dominate over the movements towards 
lower latitudes. The latter are only visible in the form of streaks 
at lower latitudes. For comparison the same quantity, but with 
m - I replaced by m - 2, has been plotted in bottom panel of 
figure 6. The similarity between the two plots indicates that the 



time and latitude dependence of the extremal values of can 
approximately be captured in this way. 

In addition to butterfly diagrams based on sunspots also 
butterfly diagrams based on observations of the radial com- 
ponent of the magnetic field are often used. In the middle 
panel of figure 6 a theoretical butterfly diagram is shown where 
Br + |B'""'|sgn(B,.) has been plotted as function of time and lat- 
itude. A comparison of the top and the middle panels of figure 
6 indicates a noticeable phase shift between the two diagrams 
which agrees roughly with the well known observation that the 
radial component of the solar magnetic field changes its sign 
near the maximum of the solar cycle (Stix, 1976). 

3. Conclusion 

In this report an attempt is described to explore the extent to 
which the operation of the solar dynamo can be understood on 
the basis of a minimal, but physically consistent, convection- 
driven dynamo model. Although the Boussinesq approxima- 
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Figure 5: (color online) Time series of magnetic dipolar energy densities and 
kinetic energy densities in the upper and lower panel respectively. The compo- 
nent Xp is shown by solid black line, while X,, Xp, and X, are shown by red, 
green and blue lines, respectively. X stands for either Mjip or E. The parameter 
values are P = 1, T = 2000, R = 140000, P„, = 3.5,/3 = 1. 

tion is a highly unrealistic assumption in the case of the Sun, 
the convection columns are similar to those found in numerical 
simulations based on anelastic models, see, for instance, Brun 
et al. (2004) and Ghizaru et al. (2010). Although there is little 
solar evidence for this type of convection, there are generally 
believed to exist as "giant cells" in the deeper region of the so- 
lar convection zone. 

The structure of the magnetic field found in our simulations 
differs substantially from the commonly assumed structure of 
the solar magnetic field dominated by a strong axisymmetric 
toroidal component. When the evidence for this traditional 
view is examined, however, it is found that observations do not 
contradict fields dominated by m - I or m - 2 components 
such as those shown in figures 2 and 5. On the contrary, the 
presence of active longitudes on the sun indicates at least that 
components with the azimuthal wavenumbers m = 1 and m - 2 
play a significant role. 

A more detailed exploration of the parameter space of our 
minimal model is likely to reveal an even better correspondence 
with solar observations. 
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Figure 6: "Butterfly" diagrams: B'^=° + \B'^=^ \sgn(B'^=°) (top), B'^'=° + 
|B^'='|sgn(B';'='') (middle), and B^=° + \B'^=^\sgniB"^=°) (bottom) are plotted 
as functions of latitude and time in the case P = 1, r = 2000, R = 140000, 
P„ = 3.5,/J= 1. 



Appendix A. Mathematical Formulation of the Numerical 
Model 

In the numerical model a spherical fluid shell rotating about 
a fixed axis described by the unit vector k is considered. It is 
assumed that a static state exists with the temperature distribu- 
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tion 

Ts ^To + ATjjr-'{l-jj)-\ 

where r denotes the distance from the center of the spherical 
shell measured in terms of multiples of the shell thickness d 
and 77 denotes the ratio of inner to outer radius of the shell. Ar 
is the temperature difference between the two boundaries. The 
gravity field is given by g = -dyr. In addition to d, the time 
d^/v, the temperature v^/yad'^, and the magnetic flux density 
v(figy^^/d are used as scales for the dimensionless description 
of the problem where v denotes the kinematic viscosity of the 
fluid, K its thermal difFusivity, g its density and fi is its magnetic 
permeability. The equations of motion for the velocity vector u, 
the heat equation for the deviation from the static temperature 
distribution, and the equation of induction for the magnetic flux 
density B are thus given by 

d,u+u 'Vu + Tkxu = -Vtt + 0r + V^h + B ■ VB, (A. 1) 



V ■ H = 0, (A.2) 
P(d,® + H ■ V0) = (/? r]r-\l - 7])-\ ■ u + 7^0, (A.3) 

V ■ B = 0, (A.4) 
V^B = P„,{d,B + u VB-B Vh), (A.5) 



where 3, denotes the partial derivative with respect to time t and 
where all terms in the equation of motion that can be written as 
gradients have been combined into Vtt. The Boussinesq approx- 
imation is assumed in that the density g is regarded as constant 
except in the gravity term where its temperature dependence, 
given by a = -{dg/dT)/g - const., is taken into account. The 
Rayleigh numbers R, the Coriolis number t, the Prandtl number 
P and the magnetic Prandtl number P,,, are defined by 

ajATd* 2Qd^ v v 

R^— , T= , P,„ = -, (A.6) 

VK V K A 

where A is the magnetic diflFusivity. Because the velocity field 
u as well as the magnetic flux density B are solenoidal vec- 
tor fields, the general representation in terms of poloidal and 
toroidal components can be used, 

H = V X (Vv X r) + Vw X r , (A.7) 
B = V X (V/j X r) + X r . (A.8) 

Equations for v and w are obtained by multiplication of the curl^ 
and of the curl of equation (1) by r. Analogously equations for 
h and g are obtained through the multiplication of equation (4) 
and of its curl by r. 

No-slip boundary conditions are used at the inner boundary 
and stress-free conditions are applied at the outer boundary, 
while the temperature is assumed to be fixed at both boundaries, 

y = 5rV = w = = O at r = r; s 77/(1 -/;) (A.9) 
V = 5^^v = 5,(w/r) = = atr = To s 1/(1 -77). (A.IO) 

Only the value 77 = 0.65 is used in the present paper. An excep- 
tion from the stress-free boundary at r = will be assumed for 
the axisymmetric part w of w, 

dr(w/r) ^ -f^/r. (A. 11) 



For the magnetic field, electrically insulating boundaries are as- 
sumed such that the poloidal function h must be matched to the 
function h^''\ which describes the potential fields outside the 
fluid shell 

g=h- = drih - /z*"') = at r = r/ and at r = r^. (A. 12) 

The numerical integration proceeds with the pseudo-spectral 
method as described by Tilgner (1999) which is based on an 
expansion of all dependent variables in spherical harmonics for 
the 9, (/3-dependences, i.e. 

V = ^ Vi'ir, t)P'l'{cos 6) expiimtp) (A. 13) 

and analogous expressions for the other variables, w, 0, h and 
g. Here P"' denotes the associated Legendre functions. For the 
r-dependence expansions in Chebychev polynomials are used. 
For the computations to be reported in this paper a minimum 
of 41 collocation points in the radial direction and spherical 
harmonics up to the order 128 have been used. 

The magnetic energy density components of dynamo solu- 
tions are defined as 

Mp = V X (VX X r) \\ M, = 1(1 Vf x r \\ 

Mp = V X (V/7 x r) M, = i(l Vg x r \\ 

where (■) indicates the average over the fluid shell and h refers 
to the axisymmetric component of h, while h is defined hyh- 
h - h. The corresponding kinetic energy densities Ep, Et, Ep 
and E, are defined analogously with v and w replacing h and g. 
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